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I describe the foundation of a Density Functional Theory approach to in- 
clude pairing correlations, which was applied to a variety of systems rang- 
ing from dilute fermions, to neutron stars and finite nuclei. Ground state 
properties as well as properties of excited states and time-dependent 
phenomena can be achieved in this manner within a formalism based on 
microscopic input. 

1. Why use a Density Functional Theory approach? 

The calculations of the ground and excited state properties of complex nu- 
clei represent an ongoing challenge for several decades. Qualitatively we 
know a lot about these nuclear properties, however, extracting these prop- 
erties from a microscopic approach is still one of the most difficult problems 
in theoretical physics. The atomic nucleus presents a number of difficulties, 
not all of them shared by other quantum many-body systems. The exis- 
tence of magic numbers is undoubtedly related to the existence of a well 
defined meanfield, in which nucleons move inside the nuclear medium. The 
experimental confirmation of the existence of superdeformed nuclei and the 
semi-quantitative explanation of the equilibrium shape deformations based 
on scmiclassical arguments represent another strong argument in favor of 
the existence of a well defined meanfield even in strongly deformed nuclei. 
Pairing correlations, which energetically represent a rather small decora- 
tion of the bulk nuclear binding energy, can equally well be described in 
terms of a generalized meanfield. At the same time the interaction among 
nucleons is strong and the description of the nuclear binding and structure 
in terms of the bare interaction among nucleons is still a formidable chal- 
lenge for the many-body theory, apart form brute force methods (quantum 
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2 A. Bulgac 

Monte Carlo (QMC) and no-core shell model) capable of handling so far 
only very light nuclei. In condensed matter theory a somewhat parallel 
approach was developed, in which the existence of well defined fermionic 
quasiparticle excitations played a fundamental role, the so called Landau 
Fermi liquid theoryjl^l^ subsequently extended to nuclei by Migdal and col- 
laboratorsP Landau took it as a given that the ground state properties 
of strongly interacting systems cannot be calculated and have to be con- 
structed in a phenomenological manner. In order to describe excited states 
one had to introduce additionally residual interactions and the emerging 
formalism was equivalent to the linear response theory.^ 

In the 1970's a serious effort was mounted to calculate nuclear properties 
on an almost theoretically self-consistent approach, using effective interac- 
tions based on a G-matrix approach to nuclear matter .13 The approach 
was cumbersome, the theoretical errors were hard to evaluate, especially 
since a number of shortcuts (based on the current dominating prevailing 
theoretical "intuition" ) were adopted. The use of the term "effective" be- 
came so widespread that one can hardly find a universal meaning, often 
being used just as a misnomer for a semi-phenomenological approach of 
one kind or another. Since the late 1970's and early 1980's we witness the 
powerful rise of the pure phenomenological approach to calculating ground 
and excited state nuclear properties using either Skyrme,'^ GognjJ^ or rela- 
tivistically inspired^ type of "effective" in-medium nuclear forces, with an 
immense proliferation of various parameterizations. There were a number 
of little nagging problems with all these models. The Skyrme interaction 
was formally interpreted as "effective" in-medium nucleon-nucleon interac- 
tion, which can be used in a Hartree-Fock like calculation.'^ Certain con- 
tributions arising from this kind of interaction lead to undesirable effects 
especially in deformed nuclei and were unceremoniously dropped. Skyrme 
parametrization of the interaction was used only in the particle-hole channel 
(similarly in relativistic models), and a totally separate phenomenological 
approach was used in the particle-particle channel, in order to treat the 
pairing correlations. Theoretically one can bring hand- waving arguments 
that such an approach was meaningful, as a similar conclusion was reached 
in the more formal Landau Fermi liquid approach,'^ where one can clearly 
show that one obtains different contributions to the irreducible diagrams 
in particle-hole and particle-particle channels respectively, even though a 
number of sum rules leads to certain correlations among them. At the same 
time Gogny and his followers insist in using the same two-body "effective" 
interaction in both the particle-hole and particle-particle channels, even 
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though it is hard to make a many-body theory argument in favor of such 
an approach. 

In the 21-st century the prevaihng theoretical attitude changed. It was 
realized and also became widely (though not unanimously) accepted that 
"effective" in-medium nuclear forces of one kind or another lack a clear 
theoretical underpinning (in particular these "effective" in-medium forces 
are not observables and moreover cannot be defined uniquely), that their 
use was merely a means to obtain an Energy Density Functional (EDF) 
(which, however, can be phenomenologically parameterized directly), and 
that an entirely new theoretical concept makes more sense, namely that of 
an Density Functional Theory (DFT), formally introduced by Kohn, Ho- 
henberg and Sham in the 1960's for electron systems-^^l DFT is in principle 
an exact approach, in which the role of the many-body wave function is 
replaced with the one-body density distribution, with the caveat that one 
needs to find an accurate energy density functional. The solution of the 
many-body Schrodinger equation is thus replaced with a significantly much 
simpler meanfield-like approach. DFT, particularly in its so called Local 
Density Approximations (LDA) of Kohn and Sham, proved to be widely 
successful in chemistry and condensed matter calculations of normal sys- 
tems, where pairing correlation are absent. Fayans^" was perhaps the first 
to introduce DFT into nuclear structure calculations in a spirit very sim- 
ilar to condensed matter physics, namely by fitting an LDA functional to 
QMC results for infinite homogeneous matter, and adding a small number 
of phenomenological gradient correction terms. 

2. DFT for a system with pairing correlations 

The DFT extension to superfluid systems has been performed in several 
ways so far, but by following a similar idea: one needs to add the anomalous 
density in order to describe the presence of a new order parameter. The 
first extension due to Gross and coUaboratorJ^i' added a dependence of 
the EDF on the nonlocal anomalous density matrix. In this way the great 
advantage of the LDA was nullified, as in this case one has to solve nonlocal 
or integro-differential equations. Especially in the case of nuclear systems, 
where pairing correlations are relatively weak and the size of the Cooper 
pair is larger than the radius of the nucleon-nucleon interaction, this kind 
of approach looks like an overkill. One expects that a large size and a 
weakly bound two-fermion Cooper pair could be accurately described using 
a zero range interaction. One can easily show however that in such a case 
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the emerging equations lead to an ultraviolet (UV) divergent anomalous 
densityP ^hat was the main reason why Gross and collaboratorJ^il 
resorted to a non-local extension of DFT to superfluid fermion systems. A 
natural UV-cutoff in electronic systems is provided by the Debye frequency. 
In nuclear systems the situation is quite different, as there are no phonon 
induced pairing correlations. 

A similar reasoning for dealing with UV-divergence was implicit as well 
in the Gogny parametrization,^ which uses a finite range "effective" in- 
medium interaction, where the finite radius of the interaction provides the 
needed UV-cutoff in the particle-particle channel. In the particle-hole chan- 
nel a finite range can be converted quite accurately (particularly in the case 
of exchange terms) into a number contact terms with spatial derivativesj^lSl 
using the density matrix expansion approach. In Skyrme-like EDFs used in 
literature most of the time practitioners introduce an arbitrary cutoff pa- 
rameter, which is often used as an additional phenomenological parameter. 
The nuclear pairing gap is small and the size of the nuclear Cooper pair 
is rather large in comparison with the nucleon-nucleon interaction radius, 
and it is hard to make the case that a finite range of the interaction is 
needed and that is responsible for the stabilization of the calculations. The 
theoretical situation here is totally similar to the nucleon-nucleon interac- 
tion for energies below the 7r-thresholcP^ where a rigorous treatment using 
contact interaction terms is easy and accurate as well to implement. 

The nature of the UV-divergence in the pairing channel was clarified a 
long time agJ^and a simple regularization scheme was later introduced.'^ 
In the case of two particles interacting with a finite range interaction the 
s-wave function for either bound states or scattering states behaves as 
oc f/|ri — outside the interaction range. In the case of scattering states 
the wave function at low energies is exp(ik-r)-|-/(fc) exp(ifc|ri— r2 |)/|ri— r2 1. 
One can show that the anomalous density matrix ^'(ri, satisfies an equa- 
tion almost identical to the Schrodinger equation for two interacting par- 
ticleJl^l and that i/(ri,r2) oc A((r)/|ri ~ r2| when |ri — r2| — > (for a 
zero-range interaction, where r = (ri + r2)/2) and exactly this is the rea- 
son why the diagonal anomalous density diverges. This divergence and 
its amplitude has physical meaning and one cannot simply hide it under 
the rug by introducing a ill-defined cutoff. Fermi devised a very simple 
approach to deal with this kind of situation, without making recourse to 
either arbitrary UV-cutoffs or to fictitious finite range effects, by introduc- 
ing a pseudo-potential.'i^ Fermi's pseudo-potential approach can be imple- 
mented in a straightforward manner in treating pairing correlations .^^^ In 
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particular, in order to describe the value of the s-wave pairing gap only 
one coupling constant is needed for both protons and neutrons, as expected 
from isospin invariancel^^ 

In a parallel approach a many-body perturbative approach to EDF, us- 
ing renormalized bare NN (and NNN) interactions was suggested by Furn- 
stahl and collaborators.^^ Only a few study cases have been considered so 
far and only in the absence of pairing correlations. 

A quite successful DFT approach for fermion systems with pairing cor- 
relations has been developed and applied to a diverse sample of physical 
systemJi^E^EsHHI and was dubbed the Superfluid Local Density Approxi- 
mations (SLDA) as a direct generalization of the LDA approach of Kohn 
and Sham. In LDA/SLDA single-particle wave functions appear explic- 
itly and often this type of approach is referred to as orbital based DFT. 
In SLDA one constructs a local EDF in terms of various densities (spin 
degrees of freedom are not shown for the sake of simplicity) 

p(r)=5^K(r)p, r(r)=^|V«„(r)p, ^(r) = ^ z.„(rK(r). (1) 

n n n 

If spin-orbit interaction is present an additional density should be added. 
The sums in the definition of these densities are all performed up to an 
UV-cutoff, in order to avoid the UV-divergence of the kinetic energy and 
anomalous densities. Even though an explicit UV-cutoff appears in the 
SLDA formulation, no dependence of the observables on this UV-cutoff 
exists once this cutoff is chosen appropriately. Both kinetic energy density 
T(r) and the anomalous density i/(r) diverge in a similar fashion and their 
contribution to EDF is handled by introducing a well chosen counter term, 
and subsequently the entire formalism becomes divergence and counter term 
bee^ One can show that a unique combination 

-^^r(r)-A(r)z.*(r) (2) 

is divergence free, where meff{r) is the effective mass, A(r) = — g^jj (r)i/(r) 
is the pairing gap and geffi'r) is the renormalized position dependent cou- 
pling constant defining the strength of the pairing correlations.'^ This 
implies that both kinetic energy and anomalous densities appear in EDF 
in this combination alone. 

Since DFT does not provide any constructive recipes for the EDF, any 
suggested approach needs validation. Unlike in the case of a phenomeno- 
logical approach in this case agreement with experiment is not the proof 
that the championed approach is correct. One needs to explicitly show that 
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the solution of the Schrodinger equation for a many-body fermion system, 
in which pairing correlations are present and the corresponding DFT in- 
carnation produce the same results for the ground state energy and ground 
state one-body density distribution. Fortunately, an extremely interesting 
system, which has mesmerized theorists across most physics subfields as 
well as experimentalists in cold atom physics, exists. This is the unitary 
Fermi gas (UFG), a system of fermions with spin- up and spin-down, inter- 
acting with a zero-range interaction and an infinite scattering length. 1^2^1211 
The UFG has properties very similar to the properties of dilute neutron 
mattei'2^t2ll as envisioned by Bertsch in IQOO.'SSl In this case it is possible to 
calculate with relatively high controlled accuracy the energy and a number 
of properties of the ground state of a large series of such systems, with 
various particle numbers with spin-up and down, both in the case of homo- 
geneous systems and systems in external confining potentials. The results 
for the ground state properties of the homogeneous state - specifically the 
ground state energy, the pairing gap and the quasiparticle excitation spec- 
trum (effective mass) - are used to build the SLDA EDF. In the case of a 
UFG the only dimensional scale in the system is the inter-particle distance 
and this fact constrains the possible EDF structure. Simple dimensional 
arguments show that apart from three dimensionless constants a, (3 and 

7 the (un-renormalized) EDF of an unpolarized (iV^ — N^) UFG has the 
simple structure: 



£sLDA[n,T, v] 



2 ^ ' 10 ^ ' ni/3(r) 



(3) 



and an additional T4a:t(r)?i(i") term, for an arbitrary external potential in 
which the system might or might not reside. As discussed above, the kinetic 
energy and the anomalous densities diverge, and a well defined renormal- 
ization procedure was devised,'^2ll2l2Sl which amounts to using these two 
densities in a unique combination in SLDA functional, see Eq. ([2]). The 
infinite matter QMC calculation^^H (where Vea;t(r) = and n(r) = const.) 
provide enough information to determine the dimensionless constants a, (3 
and 7. 

An independent series of QMC calculations of a large number of systems 
with various numbers of fermions 7V-|- and in an external harmonic trapl^ 
provide results for the ground state properties of inhomogeneous systems. 
At this point one can use the SLDA functional to predict the properties of 
these finite systems, see following table for a sample of results l^^lSfil 

The degree of agreement between the QMC results for the finite inho- 
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±\. bdllipic 




clit;ig,lt;c) iUi Uillldiy icl ililUilo 111 


a, llcLl lllUlllL. LI dp 




Eqmc 


EsLDA 


Error 


(2,1) 


4.281 ±0.004 


AA17 


3.2% 


(2,2) 


5.051 ±0.009 


5.405 


7% 


(3,2) 


7.61 ±0.01 


7.602 


0.1% 


(3,3) 


8.639 ±0.03 


8.939 


3.5% 


(4,3) 


11.362 ±0.02 


11.31 


0.49% 


(4,4) 


12.573 ±0.03 


12.63 


0.48% 


(5,5) 


16.806 ±0.04 


16.19 


3.7% 


(6,6) 


21.278 ±0.05 


21.13 


0.69% 


(7,6) 


24.787 ±0.09 


24.04 


3% 


(7,7) 


25.923 ±0.05 


25.31 


2.4% 


(8,8) 


30.876 ±0.06 


30.49 


1.2% 


(9,9) 


35.971 ±0.07 


34.87 


3.1% 


(10,10) 


41.302 ±0.08 


40.54 


1.8% 


(11,10) 


45.474 ±0.15 


43.98 


3.3% 


(11,11) 


46.889 ±0.09 


45.00 


4% 



mogeneous systems and the corresponding SLDA is a measure of the abihty 
of the DFT to describe strongly interacting superfluid fermionic systems. 
One should keep in mind that the accuracy of the QMC results is cur- 
rently at the level of 5% (the least accurate quantity being the pairing gap 
Aufg{QMC) = 0.504(24)eF and AuFoiexp.) « 0.45(5)ir_F ), and that the 
QMC results for the infinite matter and finite systems were obtained by dif- 
ferent groups, using somewhat different numerical approaches. Apart form 
this reduced sample of results presented in this table many other theoret- 
ical results (thermodynamic properties, collective states, thermodynamic 
and quantum phase transitions) ^^ ^ ^^ ^ and also extensive comparisons with 
results of many experiments are available in literature. The good quality 
of the agreement shown here, along with the theoretical arguments pre- 
sented in favor of the SLDA functional above lend strong support to the 
assertion that superfluid correlations in fermionic systems can be accurately 
described within the DFT approach. The quality of the agreement between 
the QMC results for finite systems in harmonic traps and the corresponding 
SLDA results (in particular to even-odd staggering) is even more surprising, 
as one might have expected that derivative corrections might exist, specifi- 
cally a term Vn(r)p/r7i n(r), proportional to an undefined dimensionless 
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constant. One can make the argument that such a term can be expected if 
the interaction has a finite range, and also that the introduction of such a 
term would completely destroy the agreement between the SLDA and the 
corresponding QMC results for the finite UFG systems f^^'^ll 

Very natural theoretical arguments allow us to extend the validity of 
the SLDA functional even further. By invoking local Galilean covariance 
one can show that in the un-renormalized SLDA functional ^ one has to 
add to perform the replacement 



a- — T(r) - — r(r) + (a — 1)- — 



n(r) 



(4) 



where j((r) is the one-body current density.'^iSlMI straightforward exten- 
sion exists as well to the case of a polarized UFG, when iV^ ^ A^^.E^Ml By 
changing the ratio N^/N^ a UFG undergoes a number of quantum phase 
transitions, from a uniform superfluid to a Larkin-Ovchinnikov phase, in 
which the order parameter oscillates in space, and further to a superfluid 
with relatively weak p-wave pairing and a normal state at very small tem- 
peratures.'^ The SLDA extension to include currents allows the description 
of excited states, in particular vortices, and of time-dependent phenomena. 
The existence of a DFT extension to time-dependent processes has been 
proven for quite some time.''^ The TDSLDA extension allowed so far the 
study of a large number of phenomena in a UFG system: the excitation of 
the Anderson-Higgs modes in which the magnitude of the pairing field is 
excited with a large amplitude; the Anderson-Bogoliubov sound modes; the 
vortex structur^^ni ^nd the dynamical generation of vortices and their non- 
trivial dynamics, in particular the first microscopic description of vortex 
crossing and reconnection in a fermionic superfluid leading to quantum tur- 
bulence predicted by Feynman in 1956; the generation of quantum shock 
waves and domain walls in the collision of two UFG clouds.'^ 



3. Nuclear DFT 

The nuclear EDF is slightly more complicated, as one has to embed the de- 
pendence on proton and neutron densities.'^^'^ The general principle how- 
ever is very similar, the nuclear EDF has to satisfy all required symmetries: 
rotational and translational invariance, parity and isospin symmetry, gauge 
symmetry and Galilean invariance. We will not discuss here questions re- 
lated to symmetry restoration, quantization of large amplitude motion, and 
the DFT stochastic extension, which were addressed recently elsewhere.'^ 
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An un-renormalized local nuclear EDF should have the following structure: 
£sLDA = ^;;;-'^pW + ^^^"W +£Jv[Pp(r ),Tp(r),/0„ (r),T„(r ),...] 

+ £5[Pp(r) + p„(r)](|i^p(r)p + \ua{r)\^) 

+ ss'[pp{r) + Pn{r)]{pp{r) - pjr)){\up{r)\'^ - |i^„(r)p) 

where the isospin symmetry is broken by the small difference between the 
proton Trip and neutron rn„ masses and the Coulomb interaction. The 
Coulomb interaction has a Hartree contribution in terms of the charge den- 
sity pch{^) (which is different from Pp{v) mostly due to the finite proton size) 
and an exchange-correlation contribution exc[pp{^), /On(r)], which is neither 
a simple Fock term nor its Slater approximation. Phenomenological studies 
of nuclear mass formulas show that a better fit can be typically obtained 
by neglecting the Coulomb exchange)^ and a many-body analysis of the 
Coulomb exchange in a nuclear medium show that its magnitude is signifi- 
cantly reduced^ - and this effect is encoded in the term ea;c[/Op(r), /9„(r)]. The 
term e jv [Pp (r) , (r ) , Pn (r) , t„ (r) , . . .] , where the ellipses stand for other den- 
sities such as spin densities and derivatives of pp.„(r), should be a symmetric 
function of the proton and neutron densities in order to satisfy the isospin 
invariance. The first five terms should be chosen so as to describe cor- 
rectly the QMC results for infinite neutron and symmetric nuclear matter. 
This is how Fayans has defined the first implementation of the Kohn-Sham 
DFT approach to nuclei.'^ In the case of infinite matter various derivative 
terms give a vanishing contribution and their contribution to the nuclear 
EDF has to be determined from the properties of finite nuclei. The term 
£s[Pp(r) -t- p„(r)](|fp(r)|-^ -I- |i/„(r)|2) is the most important term describing 
the nuclear pairing correlations, which also satisfies isospin invariance and 
it has been used to describe odd-even effects in more than 200 nuclei with a 
simple constant volume pairing {es = cx)nst.)W^^ The same formalism was 
used to demonstrate that a vortex in neutron matter develops a large den- 
sity depletion at its core, a feature which controls the pinning mechanism.l^Sl 
The additional term e s'[pp{v) + pn{v)]{pp{Y) - p„(r))(|i/p(r)p - |j^„(r)P) is 
also isospin symmetric, but at this time it is not entirely clear whether 
such a term is present and whether its presence is required. In many 
phenomenological nuclear mass studie^^H various authors introduce inde- 
pendent coupling constants in the particle-particle channel for neutron and 
protons, which is a clear violation of the isospin symmetry. What is even 
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more puzzling is the fact that these phenomenological studies claim that 
the strength of the coupling is stronger in the proton channel than in the 
neutron channel, which is hard to substantiate microscopically. Another 
inconsistency of these phenomenological approaches is the introduction of 
different coupling strengths for even and odd systems which is not needed 
and is theoretically unjustifiedli^ It is natural to expect that the pairing 
interaction is weaker in the proton channel, due to the repulsive character of 
the Coulomb interaction!^ Another problem with approaches which try to 
be more microscopic in this respeci"^^ "^^ is the adoption of the pairing gaps 
calculated in the weak coupling BCS approximation, which are significantly 
higher than the actual values of the pairing gaps in infinite matter!^ It is 
known since 196ll^that pairing gaps are reduced by a factor {4e)^''^ « 2.2, 
when one accounts for the contribution of the induced interactions, which 
is also confirmed by full QMC calculations of the UFG and of the neutron 
matter.'^HHH Unfortunately there are no QMC calculations of the proton and 
neutron pairing gaps in symmetric nuclear matter, and thus the possible 
dependence of these gaps on the isospin composition of the nuclear matter 
is still unknown. The divergences of the kinetic energy and anomalous den- 
sities in ^ should be dealt with as described above,^ ^^'^^ ' ^^ ' see Eq. ([2|) . 
The Galilean invariance is retrieved by using the recipe described above in 
the case of kinetic energy densitjff^HH including proton and neutron cur- 
rent densities, see Eq. ^ and similar extensions for other densities.^ With 
the inclusion of currents the DFT can be used to describe excited states and 
time-dependent phenomena and collisions .I2B |23 1 25 (33] r^j^g gj.g^ application 
of the TDSLDA extension to a nuclear process was described recentlj^^Sl by 
calculating for the first time the excitation of the giant dipole resonances in 
deformed open-shell heavy nuclei without any unjustified approximations 
and the agreement with experimental data is very good, without the need 
of any fitting parameters. Formally the solution of the TDSLDA equations 
appears as a time-dependent Hartree-Bogoliubov formalism in 3D and the 
complexity of these equations requires the use of leadership class computers, 
as it amounts to solving tens to hundreds of thousands of coupled nonlin- 
ear time-dependent 3D partial differential equations for tens to hundreds 
of thousands of time steps. With TDSLDA a microscopic treatment of low 
energy nuclear collisions, analogous to the collisions of superfluid atomic 
clouds performed recently,'^^! and induced nuclear fission in particular is 
within reach for the first time. 
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